scaleProblem <- function (T, S0, n0, y0, Proc, Y) {
  # For each year y in Y
  sapply (Y, function (y) {
    # Build the year y system
    S <- newSysPeak (Proc, peakFlops (y))

    # Determine the year y problem size
    n <- calcNewProblemSize (T, S0, n0, y0, S, y)

    return (n)
  })
}

y0 <- 0
n0 <- 4096
Y <- seq (y0, y0+10)

# Build the baseline GPU system
S0 <- newSysPeak (GPU_base, peakFlops (y0))

buildGPU <- function (y) {
  newSysPeak (GPU_base, peakFlops (y))
}

buildGPUs <- function (Y) {
  # For each year 'y' in 'Y'
  sapply (Y, function (y) { buildGPU (y) })
}

T_gpu <- function (T0, S0, n0, y0, Y, T) {
  # For each year y in Y
  sapply (Y, function (y) {
    n <- scaleProblem (T0, S0, n0, y0, GPU_base, y)
    T (buildGPU (y), n, y)
  })
}

buildCPU <- function (T0, S0, n0, y0, y) {
  # Determine a problem size appropriate for year 'y'
  n <- scaleProblem (T0, S0, n0, y0, GPU_base, y)
    
  # "Build" a CPU system that completes in the same time
  P_cpu <- matchPeak (T0, GPU_base, peakFlops (y), CPU_base, n, y)
  newSysPeak (CPU_base, P_cpu)
}

buildCPUs <- function (T0, S0, n0, y0, Y) {
  # For each year 'y' in 'Y'
  sapply (Y, function (y) { buildCPU (T0, S0, n0, y0, y) })
}

T_cpu <- function (T0, S0, n0, y0, Y, T) {
  # For each year y in Y
  sapply (Y, function (y) {
    # Determine a problem size appropriate for year 'y'
    n <- scaleProblem (T0, S0, n0, y0, GPU_base, y)
    
    # "Build" a CPU system that completes in the same time
    S_cpu <- buildCPU (T0, S0, n0, y0, y)

    # Return its execution time under measure T
    T (S_cpu, n, y)
  })
}

# Fix total time
N_total <- scaleProblem (T__FFT3D, S0, n0, y0, GPU_base, Y)
N_overlap <- scaleProblem (T_overlap__FFT3D, S0, n0, y0, GPU_base, Y)

getTime <- function (s) {
  if (s == "Total")
    T__FFT3D
  else if (s == "Intra")
    T_intra__FFT3D
  else if (s == "Inter")
    T_inter__FFT3D
  else if (s == "Overlap")
    T_overlap__FFT3D
  else
    function (...) { return (-1) }
}

getTimeCalculator <- function (s) {
  if (s == "CPU")
    T_cpu
  else if (s == "GPU")
    T_gpu
  else
    function (...) { return (-1) }
}

Match <- c ()
Eval <- c ()
Arch <- c ()
Dim <- c ()
Year <- c ()
Time <- c ()
Sockets <- c ()
Peak <- c ()
PeakRatio_gpu2cpu <- c ()
SocketsRatio_cpu2gpu <- c ()

for (m in c ("Total", "Intra", "Inter", "Overlap")) {
  T_m <- getTime (m)
  for (e in c ("Total", "Intra", "Inter", "Overlap")) {
    T_e <- getTime (e)
    for (a in c ("CPU", "GPU")) {
      T_a <- getTimeCalculator (a)
      
      Dim <- c (Dim, N)
      Year <- c (Year, Y0+Y)
      Match <- c (Match, rep (m, length (N)))
      Eval <- c (Eval, rep (e, length (N)))
      Arch <- c (Arch, rep (a, length (N)))

      Time <- c (Time, T_a (T_m, S0, n0, y0, Y, T_e))

      Peak_gpu <- peakFlops (Y)
      Sockets_gpu <- sapply (Y, function (y) { getSockets (buildGPU (y)) (y) })
      
      Peak_cpu <- sapply (Y, function (y) {
                               S_cpu <- buildCPU (T_m, S0, n0, y0, y)
                               getPeak (S_cpu) (y)
                             })
      Sockets_cpu <- sapply (Y, function (y) {
                                  S_cpu <- buildCPU (T_m, S0, n0, y0, y)
                                  getSockets (S_cpu) (y)
                                })

      if (a == "CPU") {
        Peak <- c (Peak, Peak_cpu)
        Sockets <- c (Sockets, Sockets_cpu)
      } else {
        Peak <- c (Peak, Peak_gpu)
        Sockets <- c (Sockets, Sockets_gpu)
      }

      PeakRatio_gpu2cpu <- c (PeakRatio_gpu2cpu, Peak_gpu / Peak_cpu)
      SocketsRatio_cpu2gpu <- c (SocketsRatio_cpu2gpu, Sockets_cpu / Sockets_gpu)
    } # a
  } # e
} # m

Data <- data.frame (Match=factor (Match, levels=c ("Total", "Overlap", "Inter", "Intra"))
                    , Eval=factor (Eval, levels=c ("Total", "Overlap", "Inter", "Intra"))
                    , Arch=Arch
                    , Dim=Dim
                    , Year=Year
                    , Time=Time
                    , Peak=Peak
                    , Sockets=Sockets
                    , PeakRatio_gpu2cpu
                    , SocketsRatio_cpu2gpu)

qplot (Year-Y0, Time, data=Data, geom=c ("point", "line"), colour=Arch, shape=Arch, facets=Eval ~ Match) +
  scale_x_continuous ("Years since 2010") +
  scale_colour_grey ()
if (DO_PRINT) ggsave (genRunTagPDF ("all-time"))

qplot (Year, Peak, data=Data, geom=c ("point", "line"), colour=Arch, shape=Arch, facets=. ~ Match) +
  scale_y_log10 () +
  scale_colour_grey ()
if (DO_PRINT) ggsave (genRunTagPDF ("all-peak"))

qplot (Year, Sockets, data=Data, geom=c ("point", "line"), colour=Arch, shape=Arch, facets=. ~ Match) +
  scale_y_log10 () +
  scale_colour_grey ()
if (DO_PRINT) ggsave (genRunTagPDF ("all-sockets"))

##### Plot configuration information
Data2_cpu <- subset (Data, Arch=="CPU" & Year==max (Year) & Eval=="Total")
Data2_gpu <- subset (Data, Arch=="GPU" & Year==max (Year) & Eval=="Total")
Data2 <- data.frame (Match=c (as.vector (Data2_cpu$Match), as.vector (Data2_gpu$Match))
                     , Arch=c (as.vector (Data2_cpu$Arch), as.vector (Data2_gpu$Arch))
                     , Sockets=c (Data2_cpu$Sockets-Data2_gpu$Sockets, Data2_gpu$Sockets)
                     , Peak=c (Data2_cpu$Peak, Data2_gpu$Peak-Data2_cpu$Peak)
                     , SocketsRatio_cpu2gpu=c (Data2_cpu$Sockets / Data2_gpu$Sockets)
                     )

qplot (Match, Peak/EXA, data=Data2, geom="bar", stat="identity", fill=Arch) +
  scale_x_discrete (name="Time-match Criterion") +
  scale_y_continuous (name="Exaflop/s") +
  scale_fill_grey ()
if (DO_PRINT) ggsave (genRunTagPDF ("system-peak"))

qplot (Match, SocketsRatio_cpu2gpu, data=subset (Data2, Arch=="CPU"), geom="bar") +
  scale_x_discrete (name="Time-match Criterion") +
  scale_y_continuous (name="Ratio of CPU-to-GPU Socket Count") +
  scale_fill_grey ()
if (DO_PRINT) ggsave (genRunTagPDF ("sockets-ratio"))

# eof
